* D:\E\replications\BJPS2020\rawdata\cps\drops_cps.do

* cd D:\E\replications\BJPS2020\
clear all
cap log close
log using ./rawdata/cps/drops_cps.txt, text replace
display "$S_TIME  $S_DATE"
about
ado dir

set more off
foreach f in 83 84 86 87 88 89 91 92 93 94 96 97 98 99 {
	use ./rawdata/cps/CPS_cleanmatch.dta, clear
	isid yof h_seq ppos

	set more off
	preserve
		use ./rawdata/cps/taxes/data/CPS_taxsim_tukey2.dta, clear
		*set seed 12345678
		*local n = floor((100*runiform() + 1))
		*di "`n'" // 3
		cap drop _merge
		merge m:1 yof h_seq taxunit2 using ./rawdata/cps/taxes/data/CPS_taxsim_tr`f'.dta
		cap drop _merge
		drop if ppos==.
		isid yof h_seq ppos
		bys yof h_seq (ppos): keep if _n==1
		tempfile t
		save `t'
	restore
	cap drop _merge
	merge m:1 yof h_seq using `t'
	isid yof h_seq ppos
	tab yof _merge
	drop if _merge==2

	// identifiers
	gen year=yoc
	isid pid year
	gen dataset="CPS"
	gen ccode=2

	d *_tc

	// weight
	cap drop weight
	clonevar weight=adj_pwgt // pwgt

	// HH size
	egen HHsize=rowtotal(adults kids)
	gen equivalence=sqrt(HHsize)
	label var equivalence "Equivalence scale: sqrt(HHsize)"

	// income
	// Household gross income (closest we have to DI)
	cap drop h_gi
	qui gen h_gi = htotval_tc
	qui replace h_gi = hinctot_tc if h_gi==.

	// Mkt income for files 1986, 1987
	cap drop p_mkt
	qui gen p_mkt = (i51a_tc + i51b_tc + i51c_tc) // earnings, self employment, farm
	qui replace p_mkt = p_mkt + i53b_tc + i53c_tc // add in interest and dividend income
	// Different variables after 1987 file
	qui replace p_mkt = ern_val_tc + ws_val_tc + se_val_tc + frm_val_tc if yof>1987
	qui replace p_mkt = p_mkt + int_val_tc + div_val_tc + rnt_val_tc + oi_val_tc if yof>1987

	cap drop h_mkt
	qui egen h_mkt = total(p_mkt), by(h_seq yof)

	* Taxes
	label var hfiitax_r "household federal income tax"
	label var hsiitax_r "household state income tax"
	label var hfica_r  "household OASDI/HI Payroll tax"
	egen taxes=rowtotal(hfiitax_r hsiitax_r hfica_r)
	replace taxes=. if yof==2013

	* Market income
	clonevar preY=h_mkt

	* Disposable income
	* gen postY=h_gi // this is actually gross income
	gen postY=h_gi-taxes
		// HH disposable income

	* state codes (https://cps.ipums.org/cps-action/variables/STATECENSUS#codes_section)
	label define state 11 "Maine" 12 "New Hampshire" 13 "Vermont" 14 "Massachusetts" 15 "Rhode Island" 16 "Connecticut" 21 "New York" 22 "New Jersey" 23 "Pennsylvania" 31 "Ohio" 32 "Indiana" 33 "Illinois" 34 "Michigan" 35 "Wisconsin" 41 "Minnesota" 42 "Iowa" 43 "Missouri" 44 "North Dakota" 45 "South Dakota" 46 "Nebraska" 47 "Kansas" 84 "Colorado" 85 "New Mexico" 86 "Arizona" 87 "Utah" 88 "Nevada" 91 "Washington" 92 "Oregon" 93 "California" 94 "Alaska" 95 "Hawaii" 51 "Delaware" 52 "Maryland" 53 "District of Columbia" 54 "Virginia" 55 "West Virginia" 56 "North Carolina" 57 "South Carolina" 58 "Georgia" 59 "Florida" 61 "Kentucky" 62 "Tennessee" 63 "Alabama" 64 "Mississippi" 71 "Arkansas" 72 "Louisiana" 73 "Oklahoma" 74 "Texas" 81 "Montana" 82 "Idaho" 83 "Wyoming", modify
	label val state state
	* tabstat hsiitax_r , by(state ) s(min mean max)
		
	// CPI (from C:\Dropbox\Rehm\Data\OECD\cpi\cpi_jan2016\cpi.dta)
	gen cpi=.
	replace cpi=.1227898	if year==1955
	replace cpi=.1246625	if year==1956
	replace cpi=.1288281	if year==1957
	replace cpi=.132344		if year==1958
	replace cpi=.1336816	if year==1959
	replace cpi=.1356306	if year==1960
	replace cpi=.1370828	if year==1961
	replace cpi=.1387261	if year==1962
	replace cpi=.1404459	if year==1963
	replace cpi=.1422421	if year==1964
	replace cpi=.1444969	if year==1965
	replace cpi=.1488535	if year==1966
	replace cpi=.1529809	if year==1967
	replace cpi=.159516		if year==1968
	replace cpi=.1682294	if year==1969
	replace cpi=.178051		if year==1970
	replace cpi=.1856943	if year==1971
	replace cpi=.1917707	if year==1972
	replace cpi=.2036179	if year==1973
	replace cpi=.2261275	if year==1974
	replace cpi=.2468026	if year==1975
	replace cpi=.2609809	if year==1976
	replace cpi=.2779491	if year==1977
	replace cpi=.2991593	if year==1978
	replace cpi=.3328281	if year==1979
	replace cpi=.3779237	if year==1980
	replace cpi=.416981		if year==1981
	replace cpi=.4425479	if year==1982
	replace cpi=.4567645	if year==1983
	replace cpi=.4764078	if year==1984
	replace cpi=.4932995	if year==1985
	replace cpi=.5026625	if year==1986
	replace cpi=.5210829	if year==1987
	replace cpi=.5423314	if year==1988
	replace cpi=.5685097	if year==1989
	replace cpi=.5991976	if year==1990
	replace cpi=.6245734	if year==1991
	replace cpi=.6434906	if year==1992
	replace cpi=.6624842	if year==1993
	replace cpi=.6797581	if year==1994
	replace cpi=.6988282	if year==1995
	replace cpi=.7193123	if year==1996
	replace cpi=.7361275	if year==1997
	replace cpi=.7475543	if year==1998
	replace cpi=.7639111	if year==1999
	replace cpi=.7897072	if year==2000
	replace cpi=.8120257	if year==2001
	replace cpi=.8249046	if year==2002
	replace cpi=.8436309	if year==2003
	replace cpi=.8662168	if year==2004
	replace cpi=.8956053	if year==2005
	replace cpi=.9244971	if year==2006
	replace cpi=.9508699	if year==2007
	replace cpi=.9873748	if year==2008
	replace cpi=.9838642	if year==2009
	replace cpi=1			if year==2010
	replace cpi=1.031568	if year==2011
	replace cpi=1.052915	if year==2012
	replace cpi=1.068338	if year==2013
	replace cpi=1.085669	if year==2014
	replace cpi=1.086957	if year==2015

	* Income concepts
	gen nY=postY/cpi if postY>0
	gen gY=preY/cpi if preY>0
	label var nY "HH yearly net income, not equivl."
	label var gY "HH yearly market income, not equivl."
	gen nYe=nY/equivalence
	gen gYe=gY/equivalence
	label var nYe "HH yearly net income, equivl."
	label var gYe "HH yearly market income, equivl."

	** income changes
	* top- and bottom-coding at p1 and p99
	foreach v of varlist nYe gYe { // nY gY
		bys ccode year: egen PCTTOP_`v'=wpctile(`v'), p(99) weights(weight)
		bys ccode year: egen PCTBOT_`v'=wpctile(`v'), p(1) weights(weight)
		replace `v'=PCTTOP_`v' if `v'>PCTTOP_`v' & `v'<.
		replace `v'=PCTBOT_`v' if `v'<PCTBOT_`v' & `v'<.
	}

	egen group=group(ccode year), label	
	foreach v in nYe gYe {
		xtset pid year, yearly
		egen gini_`v'_2560 = inequal(`v') if age>=25 & age<=60, by(group) weight(weight) index(gini)
		egen gini_`v'_all = inequal(`v'), by(group) weight(weight) index(gini)	
		cap drop d`v'
		gen d1_`v'=d1.`v'/abs(l1.`v') if `v'<. & l1.`v'<.
		label var d1_`v' "Growth in `v'"
		gen d`v'=d1.`v'/(abs(l1.`v'+`v')/2) if `v'<. & l1.`v'<. // arc changes
		label var d`v' "Arc change of `v'"
		gen ESI25_`v'=(d`v'<=-0.25) if d`v'<.
		gen ESI50_`v'=(d`v'<=-0.5) if d`v'<.
		gen ESI10_`v'=(d`v'<=-0.1) if d`v'<.
		gen ESI00_`v'=(d`v'<0) if d`v'<.
		gen upESI25_`v'=(d`v'>=0.25) if d`v'<.
		gen upESI50_`v'=(d`v'>=0.5) if d`v'<.
		gen upESI10_`v'=(d`v'>=0.1) if d`v'<.
		gen upESI00_`v'=(d`v'>0) if d`v'<.

		gen gESI25_`v' =(d1_`v'<=-0.25) if d1_`v'<.
		gen gESI50_`v' =(d1_`v'<=-0.5) if d1_`v'<.
		gen gESI10_`v' =(d1_`v'<=-0.1) if d1_`v'<.
		gen gESI00_`v' =(d1_`v'<0) if d1_`v'<.
		gen upgESI25_`v'=(d1_`v'>=0.25) if d1_`v'<.
		gen upgESI50_`v'=(d1_`v'>=0.5)  if d1_`v'<.
		gen upgESI10_`v'=(d1_`v'>=0.1)  if d1_`v'<.
		gen upgESI00_`v'=(d1_`v'>0)  if d1_`v'<.

		if "`v'"=="nYe" local foo "HH disposable income"
		else if "`v'"=="gYe" local foo "HH market income"
		label var gini_`v'_2560 "Gini if `foo', ages 25-60"
		label var gini_`v'_all "Gini if `foo', ages all"
		foreach j in 00 10 25 50 {
			label var ESI`j'_`v' "`j'+% arc drop, `foo'"
			label var upESI`j'_`v' "`j'+% arc gain, `foo'"
			label var gESI`j'_`v' "`j'+% drop, `foo'"
			label var upgESI`j'_`v' "`j'+% gain, `foo'"
		}
	}

	// Joint incidence of MI + DI drops
	foreach j in 00 10 25 50 {
		egen compESI`j'=group(ESI`j'_gYe ESI`j'_nYe), label
		egen compgESI`j'=group(gESI`j'_gYe gESI`j'_nYe), label	
		label var compESI`j' "group(ESI`j'_gYe ESI`j'_nYe)"
		label var compgESI`j' "group(gESI`j'_gYe gESI`j'_nYe)"
	}
	label list compESI25

	foreach j in 00 10 25 50 {
		foreach m in 1 2 3 4 {
			gen RR`j'_`m'= (compESI`j'==`m') if compESI`j'<.
			gen RR`j'_g`m'=(compgESI`j'==`m') if compgESI`j'<.
		}
		label var RR`j'_1  "MI+DI drops (ESI`j'_gYe=0 + ESI`j'_nYe==0)"
		label var RR`j'_2  "MI+DI drops (ESI`j'_gYe=0 + ESI`j'_nYe==1)"
		label var RR`j'_3  "MI+DI drops (ESI`j'_gYe=1 + ESI`j'_nYe==0)"
		label var RR`j'_4  "MI+DI drops (ESI`j'_gYe=1 + ESI`j'_nYe==1)"
		label var RR`j'_g1 "MI+DI drops (gESI`j'_gYe=0 + gESI`j'_nYe==0)"
		label var RR`j'_g2 "MI+DI drops (gESI`j'_gYe=0 + gESI`j'_nYe==1)"
		label var RR`j'_g3 "MI+DI drops (gESI`j'_gYe=1 + gESI`j'_nYe==0)"
		label var RR`j'_g4 "MI+DI drops (gESI`j'_gYe=1 + gESI`j'_nYe==1)"
	}

	foreach v in nYe gYe {
		clonevar my_`v'=`v'
		clonevar sd_`v'=`v'
		label var my_`v' "`: var label `v'' (mean)"
		label var sd_`v' "`: var label `v'' (SD)"
	}

	preserve
		keep if age>=25 & age<=60
		foreach j in 00 10 25 50 {
			foreach v in nYe gYe {
				gen d`j'_`v'=d`v' if ESI`j'_`v'==1
				label var d`j'_`v' "Median drop size if ESI`j'_`v'==1"
				gen d1_`j'_`v'=d`v' if gESI`j'_`v'==1
				label var d1_`j'_`v' "Median drop size if gESI`j'_`v'==1"
			}
		}
		foreach v of varlist _all {
			local l`v': var label `v'
		}
		collapse (sd) sd_* (mean) my_* gin*_* *ESI??_* RR??_* (median) d??_?Ye d1_??_?Ye [aw=weight], by(ccode year)
		foreach v of varlist _all {
			label var `v' "`l`v''"
		}
		gen ages="25-60"
		gen dataset="CPS"
		compress
		note: ./rawdata/cps/drops_cps.do
		note: ./rawdata/cps/draws/ESIs_ages25-60_cps_`f'.dta
		note: Extract from Stuart Craig / Austin Nichols, results from PhR
		note: Created on `= c(current_date)'
		saveold ./rawdata/cps/draws/ESIs_ages25-60_cps_`f'.dta, replace
	restore

	* collapse to get drops
	preserve
		foreach j in 00 10 25 50 {
			foreach v in nYe gYe {
				gen d`j'_`v'=d`v' if ESI`j'_`v'==1
				label var d`j'_`v' "Median drop size if ESI`j'_`v'==1"
				gen d1_`j'_`v'=d`v' if gESI`j'_`v'==1
				label var d1_`j'_`v' "Median drop size if gESI`j'_`v'==1"
			}
		}

		foreach v of varlist _all {
			local l`v': var label `v'
		}
		collapse (sd) sd_* (mean) my_* gin*_* *ESI??_* RR??_* (median) d??_?Ye d1_??_?Ye [aw=weight], by(ccode year)
		foreach v of varlist _all {
			label var `v' "`l`v''"
		}
		gen ages="all"
		gen dataset="CPS"
		compress
		note: ./rawdata/cps/drops_cps.do
		note: ./rawdata/cps/draws/ESIs_cps_`f'.dta
		note: Extract from Stuart Craig / Austin Nichols, results from PhR
		note: Created on `= c(current_date)'
		saveold ./rawdata/cps/draws/ESIs_cps_`f'.dta, replace
	restore
}

****************
* Average the draws
****************

* ages 25-60
clear
foreach f in 83 84 86 87 88 89 91 92 93 94 96 97 98 99 {
	append using ./rawdata/cps/draws/esis_ages25-60_cps_`f'.dta
}
ds ccode year ages dataset, not
foreach v of varlist _all {
	local l`v': var label `v'
}
collapse (mean) `r(varlist)', by(ccode year ages dataset)
foreach v of varlist _all {
	label var `v' "`l`v''"
}
drop if ccode==.
compress
note: Created on `= c(current_date)'
note: ./rawdata/cps/drops_cps.do
note: ./rawdata/cps/ESIs_ages25-60_cps.dta
saveold ./rawdata/cps/ESIs_ages25-60_cps.dta, replace

* all ages
clear
foreach f in 83 84 86 87 88 89 91 92 93 94 96 97 98 99 {
	append using ./rawdata/cps/draws/esis_cps_`f'.dta
}
ds ccode year ages dataset, not
foreach v of varlist _all {
	local l`v': var label `v'
}
collapse (mean) `r(varlist)', by(ccode year ages dataset)
foreach v of varlist _all {
	label var `v' "`l`v''"
}
drop if ccode==.
compress
note: Created on `= c(current_date)'
note: ./rawdata/cps/drops_cps.do
note: ./rawdata/cps/ESIs_cps.dta
saveold ./rawdata/cps/ESIs_cps.dta, replace

display "$S_TIME  $S_DATE"
cap log close


